# -*- coding: utf-8 -*-
"""
Created on Tue May  7 13:25:21 2019

@author: Solydo
"""
import numpy as np

def LogisticRegression(X,Y,coe=1E-5):
    # 梯度下降法对率回归
    # 输入X和Y，维度分别为m*n和m*1
    # coe为终止条件，当梯度小于coe时停止迭代
    # 返回拟合到的参数theta,其中最后一个参数theta[-1]为b
    m,n=np.shape(X)
    X1=np.c_[X,np.ones(m)]    #用c_函数实现数组的横向拼接
    theta=np.zeros(n+1)    #在matlab中这样生成的是 (n+1)*(n+1)的矩阵，
                           #但是在numpy中为1*(n+1)的向量
                           #注意这里向量不分行向量和列向量的
    dL=np.ones(n+1)        #梯度，初始任意设置为较大向量
    while np.linalg.norm(dL)>coe:
        p1=1/(1+np.exp(-np.dot(X1,theta)))          #计算p1,对应于(3.23)式，也就是logistic函数
        dL=-np.dot(X1.T,Y-p1)                       #梯度,结果为(n+1)维向量
        H=np.dot(np.dot(X1.T,np.diag(p1*(1-p1))),X1)   #hessian矩阵，结果为(n+1)*(n+1)
        theta-=np.dot(np.linalg.pinv(H),dL)         #更新theta
    return theta                                    #返回最佳参数

def AnalysisData(X,DataType=None):
    # 分析数据
    # DataType 返回每个属性的数据类型，若为字符串类型，说明是离散型，并统计各个可能取值
    # X1 返回数字化的数据，将离散型数值转变为相应取值的编号
    m=len(X)     #样本数
    n=len(X[0])  #属性数
    if DataType==None: #如果未从外部输入DataType,则需自动创建
        DataType=[]
        for p in range(n):
            DataType.append([type(X[0][p])])
            if DataType[p][0]==str:
                for q in range(m):
                    if X[q][p] not in DataType[p]:
                        DataType[p].append(X[q][p])
    X1=[]   #数字化的数据
    for p in range(m):
        temp=[]
        for q in range(n):
            if DataType[q][0]==str:
                temp.append(DataType[q].index(X[p][q]))
            else:
                temp.append(X[p][q])
        X1.append(temp)
    return DataType,X1

def OneHot(X,DataType,order):
    # one-hot独热编码，将原数据X变为独热编码的数据X1
    if order==None:
        order=[True]*len(DataType)
    num=[len(i)-1 for i in DataType]  #每个属性的取值数
    expand=[num[i]*(not order[i])+order[i] for i in range(len(num))] #每个属性需要扩展的数目
    X1=[]
    for i in range(len(num)):
        if expand[i]==1:
            X1=np.r_[X1,X[:,i]]
        else:
            for n in range(expand[i]):
                X1=np.r_[X1,(X[:,i]==(n+1))]
    return X1.reshape(sum(expand),-1).T

def most(Y):
    # 统计样本集中的多数类
    Yvalue=[]  #Y取值
    Ynum=[]    #取值数目
    for y in Y:
        if y not in Yvalue:
            Yvalue.append(y)
            Ynum.append(1)
        else:
            Ynum[Yvalue.index(y)]+=1
    try:
        Yvalue[Ynum.index(max(Ynum))]
    except:
        print(Yvalue,Ynum)  
    return Yvalue[Ynum.index(max(Ynum))]

def Index(Y,rule='InfoGain'):
    #计算信息熵或者基尼系数
    pass 

def score(X,Y,considerfeature,DataType,rule='InfoGain'):
    pass

def TransConsider(considerfeature,DataType,order):
    #======== 输入 ========
    # 当前考察属性(属性索引号) considerfeature
    # 各个属性的数据信息 DataType
    # 各个属性是否为有序 order
    #======== 输出 ========
    # 对无序属性的数据进行one-hot编码之后
    # 所考察属性新的索引号 newcf
    # 新的索引号所对应的原索引号 IndexBack
    num=[len(i)-1 for i in DataType]  #每个属性的取值数
    expand=[num[i]*(not order[i])+(order[i]) for i in range(len(num))] #每个属性需要扩展的数目
    newcf=[]
    IndexBack=[]
    for cf in considerfeature:
        IndexBack+=[cf]*expand[cf]
        newcf+=list(np.arange(expand[cf])+sum(expand[:cf]))
    return newcf,IndexBack

def choose(X,Y,FeatureName,considerfeature,DataType,rule='InfoGain',order=None):
    scores=[]
    if rule=='LogisticReg':
        global scale,bias
        X1=(OneHot(X,DataType,order)-np.array(bias))/np.array(scale)  #归一化
        newcf,IndexBack=TransConsider(considerfeature,DataType,order) #所考察属性的索引号需要转换
        values=list(abs(LogisticRegression(X1[:,newcf],Y)[:-1]))
        bestvalue=max(values)
        index=values.index(bestvalue)
        f=IndexBack[index]
        n=1
        while n<len(values):
            if DataType[f][0]==str:
                return [f]
            else:
                theta=LogisticRegression(X[:,[f]],Y)
                t=-theta[-1]/theta[0]
                if t<=min(X[:,[f]]) or t>=max(X[:,[f]]):
                    f=IndexBack[values.index(np.sort(values)[-n])]
                    n+=1
                else:
                    return [f,t]
        return [f,t]

# 本题不作剪枝处理，用不到，保留部分
def precut(Xt,Yt,Xv,Yv,feature,DataType,equalcut=True):
    # 预剪枝子函数
    pass

def aftercut(tree,Xt,Yt,Xv,Yv,FeatureName,LabelName,DataType,equalcut=False):
    # 后剪枝子函数,递归方法实现
    pass

scale=[]  # 用于对数据进行归一化操作
bias=[]   # X'=(X-bias)/scale
def CreatTree(Xt,Yt,Xv,Yv,FeatureName,LabelName,considerfeature='all',DataType=None,rule='InfoGain',cut=None,equalcut=True,deep=0,normal=None,order=None):
    #训练决策树
    #-------------初始化--------------
    if considerfeature=='all':   #初始时的'all'表示考虑所有特征
        considerfeature=list(range(len(FeatureName)))
    if deep==0:                  #初始时需要分析数据，将X数字化，并得到每个X的数据类型
        if order==None:
            order=[True]*len(FeatureName)
        mt=len(Yt)
        mv=len(Yv)
        DataType,X=AnalysisData(Xt+Xv,DataType)
        X=np.array(X)
        Xt=X[:mt,:]
        Xv=X[mt:,:]
        Yt=np.array(Yt)
        Yv=np.array(Yv)
        
        global scale,bias     #确定对率回归时的归一化操作参数 (X-bias)/scale ->X'
        X1=OneHot(X,DataType,order) #独热编码
        if normal==None:      #不做归一化
           bias=[0]*X1.shape[1]
           scale=[1]*X1.shape[1]
        if normal=='min-max': #min-max归一化
           bias=list(np.min(X1,axis=0))
           scale=list(np.max(X1,axis=0)-np.min(X1,axis=0))
        if normal=='z-score': #z-score归一化
           bias=list(np.mean(X1,axis=0))
           scale=list(np.std(X1,axis=0))

    #------ --考察几种特殊情况---------
    if (Yt==Yt[0]).all():      #若所有样本属于同一类
        return LabelName[Yt[0]]
    if len(considerfeature)==0:  #若待考察特征为空
        return LabelName[most(Yt)]
    if (Xt[:,considerfeature]==Xt[0,considerfeature]).all():  #若所有样本取值相同
        return LabelName[most(Yt)]
    #--------最优划分属性和预剪枝-------
    BestFeature=choose(Xt,Yt,FeatureName,considerfeature,DataType,rule,order)  #选择最佳划分属性
    bf=BestFeature[0]
    if cut=='pre' and precut(Xt,Yt,Xv,Yv,bf,DataType,equalcut):
        return LabelName[most(Yt)]
    #------------分支-----------------
    if len(BestFeature)==1:  #若为离散型属性
        nt=FeatureName[bf]+'=?'  #nt意为nodetext
        tree={nt:{}}
        considerfeature.remove(bf)
        for k in range(1,len(DataType[bf])):
            Xtk=Xt[Xt[:,bf]==k,:]
            Ytk=Yt[Xt[:,bf]==k]
            Xvk=Xv[Xv[:,bf]==k,:]
            Yvk=Yv[Xv[:,bf]==k]
            traink=[]; n=0
            if len(Ytk)==0:  #若该分支的训练集为空
                tree[nt][DataType[bf][k]]=LabelName[most(Yt)]
            else:
                tree[nt][DataType[bf][k]]=CreatTree(Xtk,Ytk,Xvk,Yvk,FeatureName,LabelName,considerfeature.copy(),DataType,rule,cut,equalcut,deep+1,normal,order)
    else: #若为连续型属性
        nt=FeatureName[bf]  #nt意为nodetext
        tree={nt:{}}
        t=BestFeature[1]  #划分点
        #----------'<=t'分支-----------
        Xt1=Xt[Xt[:,bf]<=t,:]; Yt1=Yt[Xt[:,bf]<=t]
        Xv1=Xv[Xv[:,bf]<=t,:]; Yv1=Yv[Xv[:,bf]<=t]
        #----------'>t'分支-----------
        Xt2=Xt[Xt[:,bf]>t,:]; Yt2=Yt[Xt[:,bf]>t]
        Xv2=Xv[Xv[:,bf]>t,:]; Yv2=Yv[Xv[:,bf]>t]
        
        if len(Xt1)==0 or len(Xt2)==0:   #如果某个分支为空
            print('注意，连续属性"%s"划分时，划分点在数据之外'%(nt))
            return LabelName[most(Yt)]   #不予划分，直接返回多少类标记      
        tree[nt]['<=%.3f'%t]=CreatTree(Xt1,Yt1,Xv1,Yv1,FeatureName,LabelName,considerfeature.copy(),DataType,rule,cut,equalcut,deep+1,normal,order)      
        tree[nt]['>%.3f'%t]=CreatTree(Xt2,Yt2,Xv2,Yv2,FeatureName,LabelName,considerfeature.copy(),DataType,rule,cut,equalcut,deep+1,normal,order)
    
    if cut=='after' and deep==0: #后剪枝
        tree,right=aftercut(tree,Xt,Yt,Xv,Yv,FeatureName,LabelName,DataType,equalcut)
    return treeimport numpy as np

def LogisticRegression(X,Y,coe=1E-5):
    # 梯度下降法对率回归
    # 输入X和Y，维度分别为m*n和m*1
    # coe为终止条件，当梯度小于coe时停止迭代
    # 返回拟合到的参数theta,其中最后一个参数theta[-1]为b
    m,n=np.shape(X)
    X1=np.c_[X,np.ones(m)]    #用c_函数实现数组的横向拼接
    theta=np.zeros(n+1)    #在matlab中这样生成的是 (n+1)*(n+1)的矩阵，
                           #但是在numpy中为1*(n+1)的向量
                           #注意这里向量不分行向量和列向量的
    dL=np.ones(n+1)        #梯度，初始任意设置为较大向量
    while np.linalg.norm(dL)>coe:
        p1=1/(1+np.exp(-np.dot(X1,theta)))          #计算p1,对应于(3.23)式，也就是logistic函数
        dL=-np.dot(X1.T,Y-p1)                       #梯度,结果为(n+1)维向量
        H=np.dot(np.dot(X1.T,np.diag(p1*(1-p1))),X1)   #hessian矩阵，结果为(n+1)*(n+1)
        theta-=np.dot(np.linalg.pinv(H),dL)         #更新theta
    return theta                                    #返回最佳参数

def AnalysisData(X,DataType=None):
    # 分析数据
    # DataType 返回每个属性的数据类型，若为字符串类型，说明是离散型，并统计各个可能取值
    # X1 返回数字化的数据，将离散型数值转变为相应取值的编号
    m=len(X)     #样本数
    n=len(X[0])  #属性数
    if DataType==None: #如果未从外部输入DataType,则需自动创建
        DataType=[]
        for p in range(n):
            DataType.append([type(X[0][p])])
            if DataType[p][0]==str:
                for q in range(m):
                    if X[q][p] not in DataType[p]:
                        DataType[p].append(X[q][p])
    X1=[]   #数字化的数据
    for p in range(m):
        temp=[]
        for q in range(n):
            if DataType[q][0]==str:
                temp.append(DataType[q].index(X[p][q]))
            else:
                temp.append(X[p][q])
        X1.append(temp)
    return DataType,X1

def OneHot(X,DataType,order):
    # one-hot独热编码，将原数据X变为独热编码的数据X1
    if order==None:
        order=[True]*len(DataType)
    num=[len(i)-1 for i in DataType]  #每个属性的取值数
    expand=[num[i]*(not order[i])+order[i] for i in range(len(num))] #每个属性需要扩展的数目
    X1=[]
    for i in range(len(num)):
        if expand[i]==1:
            X1=np.r_[X1,X[:,i]]
        else:
            for n in range(expand[i]):
                X1=np.r_[X1,(X[:,i]==(n+1))]
    return X1.reshape(sum(expand),-1).T

def most(Y):
    # 统计样本集中的多数类
    Yvalue=[]  #Y取值
    Ynum=[]    #取值数目
    for y in Y:
        if y not in Yvalue:
            Yvalue.append(y)
            Ynum.append(1)
        else:
            Ynum[Yvalue.index(y)]+=1
    try:
        Yvalue[Ynum.index(max(Ynum))]
    except:
        print(Yvalue,Ynum)  
    return Yvalue[Ynum.index(max(Ynum))]

def Index(Y,rule='InfoGain'):
    #计算信息熵或者基尼系数
    pass 

def score(X,Y,considerfeature,DataType,rule='InfoGain'):
    pass

def TransConsider(considerfeature,DataType,order):
    #======== 输入 ========
    # 当前考察属性(属性索引号) considerfeature
    # 各个属性的数据信息 DataType
    # 各个属性是否为有序 order
    #======== 输出 ========
    # 对无序属性的数据进行one-hot编码之后
    # 所考察属性新的索引号 newcf
    # 新的索引号所对应的原索引号 IndexBack
    num=[len(i)-1 for i in DataType]  #每个属性的取值数
    expand=[num[i]*(not order[i])+(order[i]) for i in range(len(num))] #每个属性需要扩展的数目
    newcf=[]
    IndexBack=[]
    for cf in considerfeature:
        IndexBack+=[cf]*expand[cf]
        newcf+=list(np.arange(expand[cf])+sum(expand[:cf]))
    return newcf,IndexBack

def choose(X,Y,FeatureName,considerfeature,DataType,rule='InfoGain',order=None):
    scores=[]
    if rule=='LogisticReg':
        global scale,bias
        X1=(OneHot(X,DataType,order)-np.array(bias))/np.array(scale)  #归一化
        newcf,IndexBack=TransConsider(considerfeature,DataType,order) #所考察属性的索引号需要转换
        values=list(abs(LogisticRegression(X1[:,newcf],Y)[:-1]))
        bestvalue=max(values)
        index=values.index(bestvalue)
        f=IndexBack[index]
        n=1
        while n<len(values):
            if DataType[f][0]==str:
                return [f]
            else:
                theta=LogisticRegression(X[:,[f]],Y)
                t=-theta[-1]/theta[0]
                if t<=min(X[:,[f]]) or t>=max(X[:,[f]]):
                    f=IndexBack[values.index(np.sort(values)[-n])]
                    n+=1
                else:
                    return [f,t]
        return [f,t]

# 本题不作剪枝处理，用不到，保留部分
def precut(Xt,Yt,Xv,Yv,feature,DataType,equalcut=True):
    # 预剪枝子函数
    pass

def aftercut(tree,Xt,Yt,Xv,Yv,FeatureName,LabelName,DataType,equalcut=False):
    # 后剪枝子函数,递归方法实现
    pass

scale=[]  # 用于对数据进行归一化操作
bias=[]   # X'=(X-bias)/scale
def CreatTree(Xt,Yt,Xv,Yv,FeatureName,LabelName,considerfeature='all',DataType=None,rule='InfoGain',cut=None,equalcut=True,deep=0,normal=None,order=None):
    #训练决策树
    #-------------初始化--------------
    if considerfeature=='all':   #初始时的'all'表示考虑所有特征
        considerfeature=list(range(len(FeatureName)))
    if deep==0:                  #初始时需要分析数据，将X数字化，并得到每个X的数据类型
        if order==None:
            order=[True]*len(FeatureName)
        mt=len(Yt)
        mv=len(Yv)
        DataType,X=AnalysisData(Xt+Xv,DataType)
        X=np.array(X)
        Xt=X[:mt,:]
        Xv=X[mt:,:]
        Yt=np.array(Yt)
        Yv=np.array(Yv)
        
        global scale,bias     #确定对率回归时的归一化操作参数 (X-bias)/scale ->X'
        X1=OneHot(X,DataType,order) #独热编码
        if normal==None:      #不做归一化
           bias=[0]*X1.shape[1]
           scale=[1]*X1.shape[1]
        if normal=='min-max': #min-max归一化
           bias=list(np.min(X1,axis=0))
           scale=list(np.max(X1,axis=0)-np.min(X1,axis=0))
        if normal=='z-score': #z-score归一化
           bias=list(np.mean(X1,axis=0))
           scale=list(np.std(X1,axis=0))

    #------ --考察几种特殊情况---------
    if (Yt==Yt[0]).all():      #若所有样本属于同一类
        return LabelName[Yt[0]]
    if len(considerfeature)==0:  #若待考察特征为空
        return LabelName[most(Yt)]
    if (Xt[:,considerfeature]==Xt[0,considerfeature]).all():  #若所有样本取值相同
        return LabelName[most(Yt)]
    #--------最优划分属性和预剪枝-------
    BestFeature=choose(Xt,Yt,FeatureName,considerfeature,DataType,rule,order)  #选择最佳划分属性
    bf=BestFeature[0]
    if cut=='pre' and precut(Xt,Yt,Xv,Yv,bf,DataType,equalcut):
        return LabelName[most(Yt)]
    #------------分支-----------------
    if len(BestFeature)==1:  #若为离散型属性
        nt=FeatureName[bf]+'=?'  #nt意为nodetext
        tree={nt:{}}
        considerfeature.remove(bf)
        for k in range(1,len(DataType[bf])):
            Xtk=Xt[Xt[:,bf]==k,:]
            Ytk=Yt[Xt[:,bf]==k]
            Xvk=Xv[Xv[:,bf]==k,:]
            Yvk=Yv[Xv[:,bf]==k]
            traink=[]; n=0
            if len(Ytk)==0:  #若该分支的训练集为空
                tree[nt][DataType[bf][k]]=LabelName[most(Yt)]
            else:
                tree[nt][DataType[bf][k]]=CreatTree(Xtk,Ytk,Xvk,Yvk,FeatureName,LabelName,considerfeature.copy(),DataType,rule,cut,equalcut,deep+1,normal,order)
    else: #若为连续型属性
        nt=FeatureName[bf]  #nt意为nodetext
        tree={nt:{}}
        t=BestFeature[1]  #划分点
        #----------'<=t'分支-----------
        Xt1=Xt[Xt[:,bf]<=t,:]; Yt1=Yt[Xt[:,bf]<=t]
        Xv1=Xv[Xv[:,bf]<=t,:]; Yv1=Yv[Xv[:,bf]<=t]
        #----------'>t'分支-----------
        Xt2=Xt[Xt[:,bf]>t,:]; Yt2=Yt[Xt[:,bf]>t]
        Xv2=Xv[Xv[:,bf]>t,:]; Yv2=Yv[Xv[:,bf]>t]
        
        if len(Xt1)==0 or len(Xt2)==0:   #如果某个分支为空
            print('注意，连续属性"%s"划分时，划分点在数据之外'%(nt))
            return LabelName[most(Yt)]   #不予划分，直接返回多少类标记      
        tree[nt]['<=%.3f'%t]=CreatTree(Xt1,Yt1,Xv1,Yv1,FeatureName,LabelName,considerfeature.copy(),DataType,rule,cut,equalcut,deep+1,normal,order)      
        tree[nt]['>%.3f'%t]=CreatTree(Xt2,Yt2,Xv2,Yv2,FeatureName,LabelName,considerfeature.copy(),DataType,rule,cut,equalcut,deep+1,normal,order)
    
    if cut=='after' and deep==0: #后剪枝
        tree,right=aftercut(tree,Xt,Yt,Xv,Yv,FeatureName,LabelName,DataType,equalcut)
    return tree